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DETERMINATION OF INTERSECTION OF PARTICLE FLIGHT PATH 
WITH SURFACE OF HELICAL TUBULAR VOID REGION 
by Millard L. Wohl and Carl D. Bogart 
Lewis Research Center 

SUMMARY 

A method for determining the intersections of a nuclear particle flight path and the 
wall of a helical tubular duct is presented. FORTRAN listings of computer subroutines 
which perform the numerical analysis involved are included. The subroutines are useful 
in any Monte Carlo transport analysis where complex void configurations of helical or 
pseudohelical geometry are involved. Their immediate purpose in this work is to serve 
as geometry subroutines to be used in conjunction with 05 -R, a Monte Carlo transport 
code developed at Oak Ridge National Laboratory. 


INTRODUCTION 

In many reactor and shield analyses which involve the transport of neutrons and 
gamma rays, it is necessary for high accuracy in local and emergent particle flux or 
current calculations to include the actual geometry of such system components as voids 
and coolant flow passages. An example of a coolant flow passage of complex geometry is 
a helical tubular duct. In performing a reactor -shield Monte Carlo transport analysis to 
determine the effects of such complex ducting on overall shield performance, individual 
neutron and gamma -ray histories are simulated on the digital computer. In the computa- 
tion of distances between successive particle collisions, it must be determined whether a 
void traversal occurs. If so, the points of intersection of the particle path and void sur- 
face must be determined, and the between-collisions path length in the material must be 
incremented by the in-void path length before a subsequent collision point is determined. 
It is possible to track correctly without incrementing path lengths by appropriate void- 
crossing distances. However, the method employed herein allows the determination of 
in-void fluxes, while the latter method does not. 




PHYSICAL MODEL 


The physical model is shown in figure 1. A slab of homogeneous shield material is 
penetrated by a helical duct. When the principles of descriptive geometry (ref. 1) are 
used, projections of the helical duct may be obtained, as shown in figure 2 where 
r = R ± s, sin a = k/z, and Z = R0. This is also shown by the following equations: 


or 



= 1 




- (sin^a)0 2 


(All symbols are defined in the appendix.) The coordinates of the projection of the center 
of the duct on a plane perpendicular to the axis of the cylinder about which the duct is 
wound (x,y-plane) are given by 


x = R cos u)z 

(la) 

y = R sin coz 

(lb) 


where u> = 2n/d and d is the helical pitch. Consider an unwound (inclined cylindrical) 
section of the helical tube cut by a plane normal to the axis of the cylinder about which the 
tube is wound. The elliptical intersection of the plane and inclined cylinder is shown in 
sketch (a). 
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This elliptical intersection is shown in an end view in sketch (b) . 


V 



The equation of the ellipse in the u, v-coordinate system is 

k_g>- V y (2a) 

a 2 a 2 esc 2 a 
or 

(u - R) 2 + V 2 sin 2 a = a 2 (2b) 

Consider, now, winding the inclined cylinder about an inner cylinder, as shown in fig- 
ure 2. As the tube is wound around the cylinder, vertical chords CC' of the ellipse 
(sketch (b)) transform into arcs of circles of radius r. Thus, the ” winding transforma- 
tion” is represented by 


u 


(3a) 


v ►rfl 


(3b) 


Equation (2b) now becomes 


(r - R) 2 + r 2 0 2 sin 2 o! = a 2 


(4) 
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which is the equation of the transformed elliptical, or bean-shaped, intersection shown in 
sketch (c). 


V 



where 9 is the angular displacement of points on the periphery of the ’’bean” from the 
semiminer axis of the unwrapped ellipse. Equation (4) describes the intersection of the 
helical duct surface with the x, y -plane. The shape of the intersection with any plane of 
constant z is invariant, and its position on an arbitary plane is determined by an angu- 
lar displacement of o»z. Equation (4) may thus be generalized for an arbitrary plane of 
constant z as follows: 


(r - R) 2 + R 2 sin 2 a;(0 - o>z) 2 = a 2 (5) 

Next, the intersection of the particle track and the helical duct must be considered. 
The equation of the straight line specifying the particle track may be written in para- 
metric form as follows (ref. 2): 


X = + (Xg - X^t 

(6a) 

y = Yi + (y 2 - y^t 

(6b) 

z = z 1 + (z 2 - Zj)t 

(6c) 


where 0 s t ^ 1. The starting point of a particle track is at x^, y^, z^ and the termi- 
nation point is at Xg, y 2 , z 2 . For simplicity of notation, let 

Ax = x„ - x. (7a) 
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Ay = y 2 - yj 


(7b) 


Az = - z ^ 


(7c) 


The transformation equations from x, y to r, 9 coordinates are 

2 2 2 
r = x + y 


6 = tan 


-i/y 


(8a) 

(8b) 


Substituting equations (6) modified by equations (7) into equations (8) yields 

r 2 = (x 1 + t Ax) 2 +(y 1 + t Ay) 2 


6 = tan 


-l/ y l + t Ay ^ 

iXj + t AXj 


(9a) 

(9b) 


Substitution of equations (9) into equation (5) yields 

)2 


| (Xj + t Ax) 2 + (y x + t Ay) 2 J ^ -R^ + 


R 2 sin 2 a 


tan 


-l/ y l + t Ay ' 
\Xj + t AX; 


— 1 2 


- o)(z 1 + t Az) 


= a 


(10) 


The intersections of the helical tubular surface and the straight-line particle track are 
determined by solving equation (10) for those values of t which satisfy 0 £ t ^ 1. 


NUMERICAL ANALYSIS 

Before directly attempting to solve equation (10), certain special cases maybe elimi- 
nated. Consider the cases illustrated in sketches (a), (b), and (c), where P and P' are 
the projections of initial and final collision points defining a particle track on a plane nor- 
mal to the axis of the cylindrical shell containing the duct. 
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The annular region shown in sketches (d) to (f) is the projection of the cylindrical 
shell containing the helical duct on a plane normal to the axis of the shell. In all three 
cases, no intersection of the particle track and the helical duct can occur. 

The equation of the outer circle of these annular projections is 

x 2 + y 2 = (R + a) 2 (11a) 

and the equation of the inner circle is 

x 2 + y 2 = (R - a) 2 (lib) 

When the values given in equations (6) and (7) are substituted into equation (11a), the in- 
tersection of the projected particle track and the larger circle is given by 

(x 1 + t Ax) 2 + (y 1 + t Ay) 2 = (R + a) 2 (12a) 

Performing the operations indicated in equation (12a) yields 

t 2 [(Ax) 2 + (Ay) 2 ] + 2t(x 1 Ax + y 1 Ay) + x 2 + y 2 - (R + a) 2 = 0 (12b) 


Let 


A = (Ax) 2 + (Ay) 2 

(13a) 

B = 2(Xj Ax + yj Ay) 

(13b) 

C = x 2 + y 2 - (R + a) 2 

(13 c) 
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The discriminant of equation (12b) is thus given as 


D = B - 4AC 

If D < 0, the situation is as depicted in sketch (a). For D > 0, consider 


( 14 ) 


, _ -0. 5(B + D 

1 — r 


1 / 2 ) 


and 


f _ -0. 5(B - D 1//2 ) 
2= ^ 


Let 


W = MAX (°> v 


t U p= MIN(1 ,t2) 


(15a) 


(15b) 


(16a) 

(16b) 


(17) 


If ^ow ^ ^p’ the situation is that of sketch (b). Let 

f(t) = At 2 + Bt + [x 2 + y 2 - (R - a) 2 J 

(The quantity in brackets is the equation of the inner circle. ) If f(tj 0 ) s 0 and 
f(t U p) — 0, the situation is as depicted in sketch (c). 

Consider the solution of equation (10). If the particle path length in the material is 
such that it would normally have its path terminated within the helical void, the upper 
value of t must be reset. Define g(t) as 


g(t) = 


(Xj + tAx) 2 + (y 1 + t Ay) 2 



— 

lx. + 1 Ax) 
' / 

J 1 . + . 



o)(Zj + t Az) 
2 


,R sin a 


( 18 ) 
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If g(t u p) < 0, set t = tg. The limiting solutions of equation (10) or, alternatively, 
equation (18) for g(t) = 0 are t^ QW and t 

It is necessary to determine the relative orientations of the flight path and helical 
void to see whether intersection actually occurs. For the values of t^ ow and t U p, there 
are corresponding values of z (eq. (6c)). From x(t up ), x(t low ), y(t up ), and y(t low ), two 
angles may be computed to determine whether the particle flight path penetrates the void. 
Consider the following projection of flight path and void -containing shell in the x, y-plane: 


y 



where Xj,y^ and ^ 2^2 are the coordinates of the points of intersection of the projected 
flight path and outer radius of the shell projection. Let 

a 

-l/ y l + \lD Ay 
a, = tan 1 1 "2 

\ X 1 + ^p Ax 

The angles a ^ and a 2 correspond to values of coz or to different z -coordinates of the 
intersection point of the path and void. Values of the z-coordinate in a certain specified 
range of z -values (i. e. , z (tq ow ) — z (^ U p)) correspond to real intersections and necessitate 
the solution of equation (10); others do not correspond to real intersections. 

The following values of z must be tested for a helical tube of 1^ turns: 

a 1 a 

Z- = — — + — “ (20a) 

co R sin a 



= tan 


-i/ y i + how Ay 
^ X 1 + how Ax 


(19a) 
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O'. 


Z 2 " 


Z 5 = 


Z 6 = 


CO 

R sin a 

; 3 = 

Z ^ + 2 7T 

: 4 = 

Z 2 + 277 

a 2 

a 

i 

CO 

T 

R sin a 

a 2 

a 


co R sin a 


= z,- + 27 r 


z 8 " z 6 + 2?r 


Zg = Z^ - 277 


z 10= Z 2 - 27r 


Z 11 “ z 5 " 2n 


z 12 ~ z 6 " 2?r 


(20b) 

(20c) 

(20d) 

(20e) 
(2 Of) 

(20g) 

(2 Oh) 
(20i) 
(20j) 
(20k) 
(20 1 ) 


It is not necessary to test the z -values which are not in the range 0 s z ^ d to see if 
they lie between z(t lQw ) and z(t ). 


SEQUENCING OF COMPUTER PROGRAMS AND SUBROUTINES 

When, in the sequence of computing in the controlling program (such as 05 -R), a 
collision is about to occur, the effect of the presence of the helical duct must be deter- 
mined. The controlling program must have a geometry subroutine, GEOM, which im- 
plements the mathematics discussed herein. GEOM calls for subroutine TRSOL, which 
determines the roots (values of t) of equation (18). (This equation is equivalent to equa- 

10 


I 



tion (10).) TRSOL calls for subroutine ITER, which refines these roots. When the colli- 
sion point has been determined, return is executed through TRSOL and GEOM to the main 
program. Function subprograms DES and GOFT, respectively, generate D of equa- 
tion (14) and g(t) of equation (18). LOOKZ is a utility function subprogram required by 
05-R, the main program. Figures 3 and 4 represent the logical flow of the subroutines. 
(All symbols are defined in the appendix.) 



Figure 3. - Fortran programs used in 05-R code. 
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Input quantities: 
MARK, XTWO, YTWO, 
ZTWO, XONE, YONE, 
ZONE, NMED, NREG 
ETA, ETAUSD BLZON 


Xj-XONE \ 

, Yj-VONE 
z, — ZONE 
ETAUSD - ETA j 


z 2 -ziwo 
DELX - x 2 - X, 
DELY — y 2 - Yi 
DELZ -z 2 - Zj 


ZLOWER 
► <Z2 
< ZUPPER 


A- 

(DELX) 2 + (DELY) 2 1 

B - 

2[x ] (DELX) + 

Y^DELYB 

C* 

<x x )2 + ( yi ) 2 

- (R^ 2 

D» 

D 2 - 4 AC 



1 



IN DAT -1 
R^CAPR + SMA 
R 2 -CAPR - SMA 
DA 2 = -i(RZ) 2 


' Read CAPR > 
SMA, SMD, XZ 
Y 2 , Z, ZLOWER 
ZUPPER 

v HD 


ALPHA - 
Tan _1 [SMD/ 


SMB- 
sin (ALPHA) 

H 


/ SMD> \ 
7 2( SMA) ? \ 
Cos (ALPHA) ' 


Write: You 
have made 
an angu- 
lar cut. 


Subroutine ' 
.TRSOL (x.y.z.N); 


Call TRSOL \ 

(TLOW, TUP, TR, NR) 


DT-\ 

J-1(2)NR 

XTWO -Xi 

+ [(DELX)DT]\ 

MARK --1 

LO V- 

-> DT - DT + TR(J + 1) - 

-►YTWO-/, 

+ [(DELY) DT] V- 

■> Z - ZLOWER 

6H / 

- TR(J) 

ZTWO -Zj 

+ [(DELZ)DT] / 

2021 


I Z -ZUPPER I 


/ T -(Z * zJ/DELZ 
ZTWO -Z 1 

XTWO -x-t + [(DELX)T] 

YTWO - + [(DELY)T] 

-Xj + [DELX (0.001 TUP)] 4 
yj-Yi + [DELY (0.001 TUP)] 
Zj-zj + [DELZ (0.001 TUP)] I 

\ M 


WX l»o 

1 + 


DES(TLOW) < 0 I 
and DESP<~0?r 


! c- 

(Xj) 2 + (yj) 



/D-^fD 


/ Tj-'ItB + D!/A 

^-(no) <^TLOW > TUP? 

4 _/ T 2 --i(B-D)/A 


\ TLOW -MAX(0, Ti) 
\ TUP -MIN (1, T 2 ) 


\ Z_ii 


(a! Subroutine GEOM 

Figure 4. - Flow diagrams for subroutines and functions. 


GOFT (TUP) 
> 0 ? / 



I* 

TLOW -0 
TUP -0.99 


/TUP -TUP +0.01 
Z- Z _1 +[(DELZ!TUP] 

















/ Function \ 

f Function \ 

( Subroutine \ 

iDESIX) J 

iGOFT(X) J 

l ITERIA.B.C.D ) J 


-^s. I Input quantities: 
L/nTlOW, TUP DA + ■ 


DT - (TUP 

TLOW/IOO)\ 

NR -1 

\ 

TMP-TUP 

+ DT 

T -TLOW - 

DT / 

X 

X 

t 

/ 


TLOT - \ Call \ 

T - DT W ITER (TLOT, t 
TUPT -T I* TUPT, TRtNRlJ 
n / xk) / 


XK - -XK \ 
NR -NR + 1 / 


T T < TMP ? > (No) 


+ COMMON VARIABLE 


tb) Subroutine TRSOL (TLOW, TUP, TR, NR) 


Output quanti-w R t \ 
ties: TR, NR 


0-K|> 


I Input quantities: 1 \ IymIN -lGOFTITMINH \ 

|TMI, TMA, T, XK" tmax1tma/ ^ YIWX J 1 GQFT(TMAX V 


T -TMIN + [tTMIN \ \ 

- TMAX) (YMINl/lYMAX \_ iy X 
* + YMINI] KONT -KONT + l/-* K0NT > 12 


Output [ 
quantity, T ‘ 


/TMAX -T X. 

f YMAX -|ARGM— (<J 


TMIN -T 
YMIN-|ARG| 


(c) Subroutine ITER (TMI, TMA, T, XKi. 


Input quantities: 

T[CAPR, SMN, SMI, Xi, 

* y 1( 2i, DELX, DELY, DELZ, 
OMEGA, SMB, A, B, C, Dm] + 


(Q)— <I> 


ARG-0 \ 
DZ - Zt [(DELZ)T] 1 
„ DX -Xi + [(DELX)T] 
DY- yi + [(DELY)T] 

SQR = V(DX) 2 + IDY) 2 j 


SMD = 0?XNo> 


ARG - |[JOMEGA)DZ] MODITPIE)} 1 
CO - cos (ARG) 

, SI -sin (ARG) 

PX -[lDX)CO] + [IDYJSI] 

PV -- [( DX)S I] + [IDY)CO] 

H J 


> PX = 0 and 
PY = 0 ? 


Input quantities, 
TEA. B, C] + 


(e) Function DESIT). 


Id) Function GOFT(T). 


DES- \ 

* I [AIT)] + BilTl + C/ 


Write: PX, 

PY, and zero 4 
in GIT) 


GOFT-ISQR -IAPR) 2 
_+ [ARGICAPRJSMB] 2 . 


Figure 4. - Concluded. 





















FORTRAN IV Listings of Computer Subprograms 


SIBFTC GEOMAA K EF , DD , FUL I ST , NODECK 
SUBROUTINE GEOM 

DATA TP I E / 6.2631853/ 

DATA INDAT/ 0/ 

COMMON/ GEOMC / MARK, XTWO, YTWQ, ZTWO, XOME, Y0N6, 
1 ZONE ,NMED,NREG, ETA, ETAUSD,BLZON 
COMMON/ CARL / CARR, SMA, SKD, XI, Yl, Zl, UELX, 
1DELY,DEL7,0MEGA,SMB,A,B,C,DA 
DIMENSION X(2),Y(?),Z(2),TR(10) 

EXTERNAL DES 
EXTERNAL GOFT 

FZ I T ( DDD , EEE , FFF ) = EEE + FFF + DDD 

SIGNT ( GGG , HHH , P P P ) = ( GGG-PPP ) * ( HHH-PPP ) 

X 1 =XONE 
Y 1 =YON F. 

Z 1 =ZONE 
E T AUSD= ET A 

I F ( INDAT. NE. 0 ) GO TO 4 

71 READ ( 5»?00)CAPR»SMA»SMD»XZ»YZ,ZZ,ZL0WEK,Z UPPER 
200 FORMAT ( 8F10. 5 ) 

I NDAT = 1 
R1=CAPR+SMA 
R2=CAPR-SMA 
DA=-.5*(R2 )**2 
IF ( SMD. E 0 . 0 . 0 ) GO 10 7 
50 OMEG A = TP 1 1: / SMD 
7 ALPHA=ATAN ( SMD/ ( T P I E*C A PR ) ) 
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C CHECK FOR TORUS 

IF(SMD.GT.l2.*SMA/COS (ALPHA))) GO TO 1 

51 WR I T E ( 6 , 100) 

100 FORMAT ( 30H0Y0U HAVE MADE AN ANNULAR CUT ) 
1 S M B = S I N (ALPHA) 

4 X 2 = X T W 0 
Y2 = Y T WO 
Z 2 = Z T W 0 
DEL X = X 2 — X 1 
DEL Y= Y2-Y 1 
DELZ=Z2-Z1 

I F ( SIGNT ( Z LOWER ,Z UPPER , Z2 )) 205,20 2,202 
205 CONTINUE 

A=DELX**2+DELY**2 
B=2.*( X1*DELX+Y1*DELY) 

C=X1**2+Y1**2-R1**2 
D=B*B~4.*A*C 
IF(D.GT.O.O) GO TO 5 

52 IF(A.NE.O.O) GO TO 53 

C A VER T ICAL IRAJEC I ORY 

54 IF(C.GE.O.O) GO TO 53 

55 CONTINUE 
C=X1**2+Y1**2-R2**2 
IF(C.LE.O.O) GO TO 53 

C IT MIGHT HIT 

56 T L 0 W = 0 • 0 
T U P = * 9 9 

6 TUP=TUP+.01 
Z=ZL+DELZ*TUP 

IF ( SIGNT ( Z LOWER, Z UP P ER, Z ) ) 209,202,202 
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209 CONTINUE 


I F ( GOFT ( TUP ) ) 6,2,2 
5 0=SQKT (D) 

T1=-.5*(B+D)/A 
T2=-.5*< B-D)/A 
TLOW=MAX 1 <0.0»T1) 

TUP=MIN1 ! 1 . » T2 ) 

I F ( TLOW . GE . TUP ) GO TO 53 

101 FORMAT! 1H0, 6E12.3 ) 

57 CONTINUE 
C=X1**2+Y1**2-R2**2 
DESP=DES< TUP) 

IF( ( DES(TLOW) .LE.O. ) .AND. ( DESP.LE.O. ) ) GO TO 53 

58 I P ( DESP .LE.O.) GO TO 2 

59 I F ( GOFT (TUP).GE.O. ) GO TO 60 
GO TO 2 

60 T U P = T 2 

C LIMITS ARE NOW TUP, TLOW 

C MUSI NOW SOLVE EON 17 

2 CALL TRSOL(TLOW,TUP,TR,NR) 

C NOW HAVE ALL THE ROOTS 

102 FORMAT! 15) 

IF(NR.EO.O) GO TO 53 

61 DT= 1 • 0 

00 3 J = 1 » NR » 2 

3 DT=DT+TR( J+1)-TR( J) 

XTW0=X1+DELX*DT 

YTW0=Y1*DELY*DT 

ZTW0=Z1+0ELZ*DT 
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202 MAKK=- 1 


Z=ZL0WER 

IF (ZTWU.LE.ZLUWER J G 0 TO 203 
204- Z=ZUPPER 
203 T= ( Z— Z 1 ) / L)E L Z 
Z T W 0 = Z 

XTW0=X1+DELX*T 
YTWU=Y1+DELY*T 
X1=X1+DELX*{ TUP+1 .E-3 ) 
Y1=Y1+DELY*( TUP+1. E-3) 
Z1=Z1+DELZ*< TUP+1 .E-3) 

GO TO 4 
33 RETURN 
END 
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$ I B F T C TRSOLS REF , DD , FUL I ST , NODECK 

SUBROUTINE TR SOL ( TLOW , T JP , TR » NR ) 

C UTILITY SR USED BY BOGART/WOHL VERS I ON' OF OBR SR# GEOM (GEOMAA) 
DIMENSION T R ( 10) 

COMMON/CARL/ CAPR»SMA»SMDtXl ,Y1 ,Z1 »DELX, 

1 DEL Y » DELZ » OMEGA » SMB t A » B » C » DA 
EXTERNAL DES 
EXTERNAL GOFT 
DT= ( TUP -T LOW ) / 100 • 

NR= I 

TMP=TUP+DT 
T=TLOW-DT 
XK = 1 . 

1 CONTINUE 

C CHECK IF NEAR ORIGIN 

IF(DES(T).LE.DA) GO TO 2 

3 I F ( (GOFT (T)*XK) .GT. 0.0) GO TO 2 

C MUST HAVE PASSED A ROOT 

4 TLOT = T — DT 
TUP T = T 

CALL ITER(TLOT t TUPT»TR(NR) »XK) 

XK=— XK 
NR=NR+1 

2 T=T+DT 

IF(T.LE.TMP) GO TO 1 

5 NR=NR- 1 
RETURN 
END 
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SIBFTC ITEKA K E F , UD , F UL I ST » NODECK 

SUBROUTINE I TER ( TM1 , TMA , T , XK ) 

C UTILITY SR USED BY BOGART/WOHL VERSION OF 05R SR# GEOM (GEOMAA) 
EXTERNAL GOFT 
K 0 N T = 0 
TMIN=TMI 
TMAX = TMA 

C METHOD OF REGULA FALSI 

YM I N = AB S (GUFT(TMIN)) 

YMAX=ABS ( GOFT ( TMAX ) ) 

A T=TMIN+(TMAX-TMIN)*YMIN/( YMAX + YMIN) 

KONT=KUNT+I 

IFIK0NT.GT.12) GO TO 2 
6 ARG = GOF TIT) *XK 

I F ( A B S (ARG).LT.l. E-6 ) GO TO 2 
5 IF(ARG)1»2,3 
3 T M I N = T 

YM1N=A8S (ARG) 

GO TO A 

1 T M A X = T 

Y MAX= AB S (ARG) 

GO TO A 

2 CONTINUE 

103 FORMAT ( 1H0» I5.3E12.3) 

RETURN 

END 
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$ I BF T C GOFT KEF , L)D, FUL I ST , NODECK 
FUNCTION GOFT(T) 

C UTILITY FUNCTION USED BY BOGART/WOHL VERSION OF 05R SK+ GEOM ( GEOMAA ) 
COMMON /CARL /CAPR,SMA,SMD»X1 , Y1 , Z1 ,DELX, 

1DELY,DELZ» OMEGA »SMB»A»B»C»DA 
DATA TP I E / 6.2831853/ 

ARG=0. 

DZ = Z 1 +DE L Z*T 
UX=X1+DELX* r 
DY=Y 1 +DEL Y*T 
SQR=SQRT ( CX**2+DY**2 ) 

C SMD=0 IMPLIES OMEG A= I NF I N I T Y 

IF(SMD.EO.O.O) GO TO 1 
2 ARG = AMOD ( OMEGA*DZ , TP I E ) 

CO=COS I ARG ) 

S I = S I N (ARG) 

PX=DX*CO+DY*SI 

PY=-DX*SI+DY*CO 

IF ( ( PY .EO. O.O) .AND. ( PX .EQ. 0.0 ) ) GO TU 99 

ARG=ATAN2 ( PY, PX ) 

GO TO 1 
99 ARG = 0.0 

WRITEl 6, 555) 

555 FORMAT < 21H PX, PY , ZERO IN G ( T ) ) 

1 GOFT = (SQR-CAPK) **2+ ( ARG*C APK*SMB ) **2 -SM A*SMA 
RETURN 
END 
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SIBFTC DES REF»DD»FULIST» NOD EC K 

FUNCTION UES(T) 

C UTILITY FUNCTION USED BY BOGART/WOHL VERSION OF 05R SR+ GEUM (GEOMAA) 
COMMON/ CARL/ CAPR, SMAtSMD, XI ,Y1 T ZI,DELX, 
lDELY,DELZ,OMEGA t SMB»A,B,C»DA 
DES=( A*T+B)*T+C 
RETURN 
END 

$ I 8F TC LOOKZ R EF , DD , FUL I ST , NODECK 
SUBROUTINE LOOKZI X» Y, Z ) 

C UTILITY FUNCTION USED BY BOGART/WOHL VERSION OF 05R SR+ GEOM (GEOMAA) 
NMED = 1 
NREG = 1 
BL Z0N = 0 . 0 
RETURN 
END 

MBFTC JOMIN REF, DDrFULIST, NODECK 
SUBROUTINE JOMINI NADD ) 

C UTILITY FUNCTION USED BY BOGART/WOHL VERSION OF 05R SR+ GEUM (GEOMAA) 
CUMMON/GEOM/HARK , X TWO , YTWO * ZT WO , XONE t YONE r 
1 ZONE t NMED, NRFG, ETA , ETAUSD » RLZON 
N ADD=0 
RETURN 
END 


Lewis Research Center, 

National Aeronautics and Space Administration, 
Cleveland, Ohio, October 19, 1967, 
126-15-01-03-22. 
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APPENDIX - SYMBOLS 


Mathematical 

FORTRAN 

Description 

symbol 

variable 



a 

SMA 

Internal radius of ducts, cm 

D 

DES 

Discriminant of eq. (12b) 

d 

SMD 

Pitch of helical duct, cm 

g(t) 

GOFT 

g(t) 

k 


See fig. 2 

l 


See fig. 2 

R 

CAPR 

Radius of cylinder around which in- 
ternal duct is wound, cm 

r 


See fig. 2 

s 


See fig. 2 

t 


Particle track parameter 

^low 

TLOW 

Lower particle track parameter 
(inner circle) 

t up 

TUP 

Upper particle track parameter 
(outer circle, sketch (b)) 

x i’ y r z i 


Coordinates of starting point of par 
tide path 

x 2’ y 2’ Z 2 


Coordinates of termination point of 
particle path 

Ax 

DELX 

Ax 

Ay 

DELY 

Ay 

Az 

DELZ 

Az 

a 

ALPHA 

Pitch angle of helical duct, deg 

J8 


See fig. 2 

9 


See fig. 2 
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OMEGA 


a> 


Z LOWER 


Z UPPER 


Parameter specifying position of in 
tersection of helical tube and 
x, y-plane 

Lower -boundary z -coordinate of 
slab of material penetrated by 
void with axis of revolution in 
z-direction 

Upper -boundary z -coordinate of 
slab 
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